Normative Jump Data Malaysia

library(readr)
library(dplyr)
library(ggplot2)
library(PupillometryR)
library(ggprism)
library(flextable)
library(officer)
library(table1)
library(janitor) 
library(caret)
library(ggpubr)
library(plotly)
library(kableExtra)
library(htmlwidgets)
library(webshot)
library(ggcorrplot)
library(factoextra)

Df <- read_csv("Project-Session-09_13_22-Session-09_13_22-AveragesAndTests_Countermovement_Jump.csv")

## Rows: 575 Columns: 88
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr  (7): TestId, Date, Time, Name, Type, Time to Stabilization, L|R Landing...
## dbl (77): System Weight, Jump Height, Jump Momentum, Countermovement Depth, ...
## lgl  (4): Segment, Position, Excluded, Tags
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

## Reduce jump trials to average
Df_average <- Df %>% filter(TestId == "AVERAGE")
Df_average <- clean_names(Df_average)
Df_average <- subset(Df_average, select = -c(test_id, date, time, segment, position, excluded, type, tags, rsi, flight_time,takeoff_velocity))
Df_average$'l_r_landing_impulse_index'<- as.numeric(Df_average$'l_r_landing_impulse_index')
Df_average$'time_to_stabilization'<- as.numeric(Df_average$'time_to_stabilization')

Data Wrangling and Cleaning

Names of participants and sex needs to be join into the main jump dataset to have a complete dataset. New dataset is then saved into an csv document

## import participant details
library(readxl)
Df_details <- read_excel("Athletes details SUKMA XX.xlsx")

# change capital name letters so that data sets match
Df_details <- Df_details %>% rename(name = NAME)

# re name males and females
Df_details$Sex <- recode_factor(Df_details$Sex,  F = "Female",
                                M = "Male")
# merge data sets by Name
Df <- inner_join(x=Df_details,y=Df_average,
                 by = "name")

## Subset NA
Df <- Df |>
  subset(!is.na(jump_height))

# convert spaces to _
Df <- clean_names(Df)

#remove names
Df <- subset(Df, select = -c(name))


write.csv(Df, "C:\\Users\\Samuel\\OneDrive\\Documents\\R\\SUKMA_Wushu\\SUKMA_Wushu_Data.csv", row.names = FALSE)

Zscore

Creation and distribuition of the zscore. First need to create two data frames, one for males and one for females.

Tscore for CMJ

next step is to create a tscore based on the zscore * 10 + 50 as recommended by the article. Then, a rank class is created based on the tscore. Lastly, factors are reorderd so that they appear in that order in the plots.

# rank based on z score
Df_male <- Df_male%>% mutate(Tscore =  zscore*10+50)
Df_female <- Df_female%>% mutate(Tscore =  zscore*10+50)

Df_male <- Df_male %>% mutate(Description = 
                            case_when(zscore > 2.0 ~ 'Excellent',
                             zscore <= 2.0 & zscore > 1.5  ~ 'Very Good',
                            zscore <= 1.5 & zscore > 1  ~ 'Good',
                            zscore <= 1 & zscore > 0.5  ~ 'Above Average',
                             zscore <= 0.5 & zscore > -0.5  ~ 'Average',
                             zscore <= -0.5 & zscore > -1.0  ~ 'Below Average',
                              zscore <= -1.0 & zscore > -1.5  ~ 'Poor',
                              zscore <= -1.5 & zscore > -2.0  ~ 'Very Poor',
                              zscore <= -2.0  ~ 'Extremely Poor'))

Df_male <- Df_male %>% mutate(Description_T = 
                            case_when(Tscore > 80 ~ 'Excellent',
                                    Tscore < 81 & Tscore > 70  ~ 'Very Good',
                                    Tscore < 71 & Tscore > 60  ~ 'Good',
                                    Tscore < 61 & Tscore > 55  ~ 'Above Average',
                                    Tscore < 56 & Tscore > 45  ~ 'Average',
                                    Tscore < 46 & Tscore > 40  ~ 'Below Average',
                                    Tscore < 41 & Tscore > 30  ~ 'Poor',
                                    Tscore < 31 & Tscore > 20  ~ 'Very Poor',
                                    Tscore < 21  ~ 'Extremely Poor')) 

Df_male$Description <- factor(Df_male$Description,
                             levels = c("Excellent", "Very Good", "Good", 
                                        "Above Average", "Average",
                                        "Below Average","Poor", 
                                        "Very Poor", "Extremely Poor"),
                             ordered = TRUE)

Df_female <- Df_female %>% mutate(Description = 
                            case_when(zscore > 2.0 ~ 'Excellent',
                             zscore <= 2.0 & zscore > 1.5  ~ 'Very Good',
                            zscore <= 1.5 & zscore > 1  ~ 'Good',
                            zscore <= 1 & zscore > 0.5  ~ 'Above Average',
                             zscore <= 0.5 & zscore > -0.5  ~ 'Average',
                             zscore <= -0.5 & zscore > -1.0  ~ 'Below Average',
                              zscore <= -1.0 & zscore > -1.5  ~ 'Poor',
                              zscore <= -1.5 & zscore > -2.0  ~ 'Very Poor',
                              zscore <= -2.0  ~ 'Extremely Poor'))

Df_female <- Df_female %>% mutate(Description_T = 
                            case_when(Tscore > 80 ~ 'Excellent',
                              Tscore < 81 & Tscore > 70  ~ 'Very Good',
                              Tscore < 71 & Tscore > 60  ~ 'Good',
                              Tscore < 61 & Tscore > 55  ~ 'Above Average',
                              Tscore < 56 & Tscore > 45  ~ 'Average',
                              Tscore < 46 & Tscore > 40  ~ 'Below Average',
                              Tscore < 41 & Tscore > 30  ~ 'Poor',
                              Tscore < 31 & Tscore > 20  ~ 'Very Poor',
                              Tscore < 21  ~ 'Extremely Poor')) 

Df_female$Description <- factor(Df_female$Description,
                             levels = c("Excellent", "Very Good", "Good", 
                                        "Above Average", "Average",
                                        "Below Average","Poor", 
                                        "Very Poor", "Extremely Poor"),
                             ordered = TRUE)

##Interactive Plots

p1 <- Df_male %>%
  ggplot(aes(y = jump_height, x = zscore, fill = Description)) +
  geom_point(position = position_jitter(width = .25), size = 4, shape = 21, color = "black") +
  scale_fill_manual(values = c("#006600", "#33CC33", "#00CC00",
                               "#ccffcc", "#FFFFFF", "#FFCC00",
                               "#CC3232", "#a32828", "#7a1e1e")) +
  theme_bw() +
  labs(title = "Male", y = "Jump Height (m)")

ggplotly(p1)

Figure 1. Interactive plot of normative jump data for male athletes.

#ggsave("p1.png")

summary_male_table <- Df_male %>%
  group_by(Description) %>%
  summarise(`Number of Athletes` = n(),
            `Min Jump Height (cm)` = format(min(jump_height), digits = 2),
            `Max Jump Height (cm)` = format(max(jump_height), digits = 2))

summary_male_table$`Jump Height (cm) Range` <- ifelse(
  summary_male_table$Description == "Excellent",
  paste("> ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Excellent"]),
  
  ifelse(summary_male_table$Description == "Very Good",
    paste(summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Very Good"],  " - ",summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Excellent"]
    ),
    ifelse(summary_male_table$Description == "Good",
            paste(summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Good"], " - ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Very Good"]),
           
        ifelse(
          summary_male_table$Description == "Above Average",
          paste(summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Above Average"], " - ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Good"]),
        
        ifelse(
          summary_male_table$Description == "Average",
          paste(summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Average"], " - ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Above Average"]),
          
               ifelse(
          summary_male_table$Description == "Below Average",
          paste(summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Below Average"], " - ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Average"]),
          
      ifelse(
        summary_male_table$Description == "Poor",
        paste(summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Poor"], " - ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Below Average"]),
        
      ifelse(
        summary_male_table$Description == "Very Poor",
        paste(summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Very Poor"], " - ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Poor"]),
            
    ifelse(
      summary_male_table$Description == "Extremely Poor",
      paste("< ", summary_male_table$`Min Jump Height (cm)`[summary_male_table$Description == "Very Poor"]),
      
          

              paste(summary_male_table$`Min Jump Height (cm)`, " - ", summary_male_table$`Max Jump Height (cm)`)
            )
          )
        )
      )
    )
  )
)
))

summary_male_table$Description <- factor(summary_male_table$Description,
                                         levels = c("Excellent", "Very Good", "Good", "Above Average", "Average", "Below Average", "Poor", "Very Poor","Extremely Poor"))

summary_male_table <- summary_male_table %>%
  select(-`Min Jump Height (cm)`, -`Max Jump Height (cm)`) %>%
  kable("html", align = "c") %>%
  column_spec(1:3, bold = TRUE) %>%
  row_spec(1:8, bold = TRUE) %>%
  kable_styling(full_width = TRUE) %>%
  row_spec(which(summary_male_table$Description == "Excellent"), background = "#00AA00") %>%
  row_spec(which(summary_male_table$Description == "Very Good"), background = "#00CC00") %>%
  row_spec(which(summary_male_table$Description == "Good"), background = "#00EE00") %>%
  row_spec(which(summary_male_table$Description == "Above Average"), background = "#99FF99") %>%
  row_spec(which(summary_male_table$Description == "Average"), background = "#FFFFFF") %>%
  row_spec(which(summary_male_table$Description == "Below Average"), background = "#FFF000") %>%
  row_spec(which(summary_male_table$Description == "Poor"), background = "#FF0000") %>%
  row_spec(which(summary_male_table$Description == "Very Poor"), background = "#CC0000") %>%
  row_spec(which(summary_male_table$Description == "Extremely Poor"), background = "#990000")

summary_male_table

Description	Number of Athletes	Jump Height (cm) Range
Excellent	2	> 0.62
Very Good	1	0.58 - 0.62
Good	9	0.55 - 0.58
Above Average	6	0.52 - 0.55
Average	32	0.47 - 0.52
Below Average	7	0.44 - 0.47
Poor	5	0.41 - 0.44
Very Poor	4	0.38 - 0.41
Extremely Poor	1	< 0.38

p2 <- Df_female %>%
  ggplot(aes(y = jump_height, x = zscore, fill = Description)) +
  geom_point(position = position_jitter(width = .25), size = 4, shape = 21, color = "black") +
  scale_fill_manual(values = c("#006600", "#33CC33", "#00CC00",
                               "#ccffcc", "#FFFFFF", "#FFCC00",
                               "#CC3232", "#a32828", "#7a1e1e")) +
  theme_bw() +
  labs(title = "Female", y = "Jump Height (m)")

ggplotly(p2)

Figure 1. Interactive plot of normative jump data for female athletes.

#ggsave("p2.png")

Notice that we don’t have a “very good” due to the zscore and the low sample size, those that section has been silenced (#), this should be improved with a larger sample size. For a larger sample remove the #, and add a closure “)”

summary_female_table <- Df_female %>%
  group_by(Description) %>%
  summarise(
    `Number of Athletes` = n(),
    `Min Jump Height (cm)` = format(min(jump_height), digits = 2),
    `Max Jump Height (cm)` = format(max(jump_height), digits = 2)
  )

summary_female_table$`Jump Height (cm) Range` <- ifelse(
  summary_female_table$Description == "Excellent",
  paste("> ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Excellent"]),
  
  #ifelse(summary_female_table$Description == "Very Good",
    #paste(summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Very Good"],  " - ",summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Excellent"]
    #),
    ifelse(summary_female_table$Description == "Good",
            paste(summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Good"], " - ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Excellent"]),
           
        ifelse(
          summary_female_table$Description == "Above Average",
          paste(summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Above Average"], " - ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Good"]),
        
        ifelse(
          summary_female_table$Description == "Average",
          paste(summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Average"], " - ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Above Average"]),
          
               ifelse(
          summary_female_table$Description == "Below Average",
          paste(summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Below Average"], " - ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Average"]),
          
      ifelse(
        summary_female_table$Description == "Poor",
        paste(summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Poor"], " - ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Below Average"]),
        
      ifelse(
        summary_female_table$Description == "Very Poor",
        paste(summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Very Poor"], " - ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Poor"]),
            
    ifelse(
      summary_female_table$Description == "Extremely Poor",
      paste("< ", summary_female_table$`Min Jump Height (cm)`[summary_female_table$Description == "Very Poor"]),

              paste(summary_female_table$`Min Jump Height (cm)`, " - ", summary_female_table$`Max Jump Height (cm)`)
            )
          )
        )
      )
    )
  )
)
)

summary_female_table$Description <- factor(
  summary_female_table$Description,
  levels = c(
    "Excellent",
    "Very Good",
    "Good",
    "Above Average",
    "Average",
    "Below Average",
    "Poor",
    "Very Poor",
    "Extremely Poor"
  )
)

summary_female_table <- summary_female_table %>%
  select(-`Min Jump Height (cm)`, -`Max Jump Height (cm)`) %>%
  kable("html", align = "c") %>%
  column_spec(1:3, bold = TRUE) %>%
  row_spec(1:8, bold = TRUE) %>%
  kable_styling(full_width = TRUE) %>%
  row_spec(which(summary_female_table$Description == "Excellent"), background = "#00AA00") %>%
  row_spec(which(summary_female_table$Description == "Very Good"), background = "#00CC00") %>%
  row_spec(which(summary_female_table$Description == "Good"), background = "#00EE00") %>%
  row_spec(which(summary_female_table$Description == "Above Average"), background = "#99FF99") %>%
  row_spec(which(summary_female_table$Description == "Average"), background = "#FFFFFF") %>%
  row_spec(which(summary_female_table$Description == "Below Average"), background = "#FFF000") %>%
  row_spec(which(summary_female_table$Description == "Poor"), background = "#FF0000") %>%
  row_spec(which(summary_female_table$Description == "Very Poor"), background = "#CC0000") %>%
  row_spec(which(summary_female_table$Description == "Extremely Poor"), background = "#990000")

summary_female_table

Description	Number of Athletes	Jump Height (cm) Range
Excellent	1	> 0.49
Good	7	0.42 - 0.49
Above Average	8	0.4 - 0.42
Average	18	0.35 - 0.4
Below Average	4	0.33 - 0.35
Poor	6	0.3 - 0.33
Very Poor	3	0.29 - 0.3
Extremely Poor	1	< 0.29

ggarrange(p1,p2)+ theme_prism()

ggsave("normative_bysex.png")

## Saving 10 x 5 in image

Distribuition of Jump Height by groups

p3 <- Df %>% ggplot(aes(x = sex, y = jump_height, fill = sex)) +
  geom_flat_violin(aes(fill = sex),
      position = position_nudge(x =.1, y=0), adjust = 1.5,trim = F, alpha =.5) +
  geom_point(aes(x=sex, y=jump_height,colour=sex),
             position = position_jitter(width = .05), size = 2.5, shape = 20)+
  geom_boxplot(aes(x = sex, y = jump_height, fill = sex),outlier.shape
               = NA, alpha = .5, width = .15, colour = "black")+
  coord_cartesian(ylim=c(0.22,0.72)) +
  scale_colour_brewer(palette = "Dark2")+
  scale_fill_brewer(palette = "Dark2")+ theme_bw() +
  ylab('Vertical Jump Height (m)')+
  xlab('Group') + theme_prism()
p3 

## Warning: Using the `size` aesthetic with geom_polygon was deprecated in ggplot2 3.4.0.
## ℹ Please use the `linewidth` aesthetic instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

p4 <- Df %>% ggplot(aes(x = events, y = jump_height, fill = sex)) +
  geom_flat_violin(aes(fill = sex),
                   position = position_nudge(x =.1, y=0), adjust = 1.5,trim = F, alpha =.5) +
  geom_point(aes(x=events, y=jump_height,colour=sex),
             position = position_jitter(width = .05), size = 2.5, shape = 20)+
  geom_boxplot(aes(x = events, y = jump_height, fill = sex),outlier.shape
               = NA, alpha = .5, width = .15, colour = "black")+
  coord_cartesian(ylim=c(0.22,0.75)) +
  scale_colour_brewer(palette = "Dark2")+
  scale_fill_brewer(palette = "Dark2")+ theme_bw() +
  ylab('Vertical Jump Height (m)')+
  xlab('Group') + theme_prism()
p4

p5 <- Df %>% ggplot(aes(x = events, y = jump_height, fill = events)) +
  geom_flat_violin(aes(fill = events),
                   position = position_nudge(x =.1, y=0), adjust = 1.5,trim = F, alpha =.5) +
  geom_point(aes(x=events, y=jump_height,colour=events),
             position = position_jitter(width = .05), size = 2.5, shape = 20)+
  geom_boxplot(aes(x = events, y = jump_height, fill = events),outlier.shape
               = NA, alpha = .5, width = .15, colour = "black")+
  coord_cartesian(ylim=c(0.22,0.8)) +
  scale_colour_brewer(palette = "Dark2")+
  scale_fill_brewer(palette = "Dark2")+ theme_bw() +
  ylab('Vertical Jump Height (m)')+
  xlab('Group') + theme_prism()
p5

ggarrange(p5,p4, ncol = 1)

table1(~jump_height + jump_momentum + 
         m_rsi + stiffness + positive_net_impulse | sex, data=Df )

	Female (N=48)	Male (N=67)	Overall (N=115)
jump_height
Mean (SD)	0.371 (0.0480)	0.494 (0.0582)	0.443 (0.0812)
Median [Min, Max]	0.376 [0.268, 0.489]	0.492 [0.358, 0.679]	0.442 [0.268, 0.679]
jump_momentum
Mean (SD)	136 (16.8)	187 (24.6)	166 (33.5)
Median [Min, Max]	137 [103, 177]	183 [127, 232]	164 [103, 232]
m_rsi
Mean (SD)	0.437 (0.0862)	0.557 (0.102)	0.507 (0.112)
Median [Min, Max]	0.437 [0.261, 0.621]	0.564 [0.245, 0.816]	0.506 [0.245, 0.816]
stiffness
Mean (SD)	-3470 (789)	-3930 (826)	-3740 (839)
Median [Min, Max]	-3300 [-5050, -1780]	-3890 [-5830, -1760]	-3670 [-5830, -1760]
positive_net_impulse
Mean (SD)	208 (28.6)	277 (39.1)	248 (48.7)
Median [Min, Max]	209 [151, 284]	272 [179, 356]	244 [151, 356]

Zscore for RSI

Creation and distribution of the zscore. First need to create two data frames, one for males and one for females.

Tscore for RSI

next step is to create a tscore for mrsi based on the zscore * 10 + 50 as recommended by the article. Then, a rank class is created based on the tscore. Lastly, factors are reorderd so that they appear in that order in the plots.

# rank based on T score
Df_male_mrsi <- Df_male_mrsi%>% mutate(Tscore =  zscore*10+50)
Df_female_mrsi <- Df_female_mrsi%>% mutate(Tscore =  zscore*10+50)

Df_male_mrsi <- Df_male_mrsi %>% mutate(Description = 
                            case_when(Tscore > 80 ~ 'Excellent',
                                    Tscore < 81 & Tscore > 70  ~ 'Very Good',
                                    Tscore < 71 & Tscore > 60  ~ 'Good',
                                    Tscore < 61 & Tscore > 55  ~ 'Above Average',
                                    Tscore < 56 & Tscore > 45  ~ 'Average',
                                    Tscore < 46 & Tscore > 40  ~ 'Below Average',
                                    Tscore < 41 & Tscore > 30  ~ 'Poor',
                                    Tscore < 31 & Tscore > 20  ~ 'Very Poor',
                                    Tscore < 21  ~ 'Extremely Poor')) 

Df_male_mrsi$Description <- factor(Df_male_mrsi$Description,
                             levels = c("Excellent", "Very Good", "Good", 
                                        "Above Average", "Average",
                                        "Below Average","Poor", 
                                        "Very Poor", "Extremely Poor"),
                             ordered = TRUE)


Df_female_mrsi <- Df_female_mrsi %>% mutate(Description = 
                            case_when(Tscore > 80 ~ 'Excellent',
                                    Tscore < 81 & Tscore > 70  ~ 'Very Good',
                                    Tscore < 71 & Tscore > 60  ~ 'Good',
                                    Tscore < 61 & Tscore > 55  ~ 'Above Average',
                                    Tscore < 56 & Tscore > 45  ~ 'Average',
                                    Tscore < 46 & Tscore > 40  ~ 'Below Average',
                                    Tscore < 41 & Tscore > 30  ~ 'Poor',
                                    Tscore < 31 & Tscore > 20  ~ 'Very Poor',
                                    Tscore < 21  ~ 'Extremely Poor')) 

Df_female_mrsi$Description <- factor(Df_female_mrsi$Description,
                             levels = c("Excellent", "Very Good", "Good", 
                                        "Above Average", "Average",
                                        "Below Average","Poor", 
                                        "Very Poor", "Extremely Poor"),
                             ordered = TRUE)

##Interactive Plots

p6 <- Df_male_mrsi %>% ggplot(y=m_rsi,x=zscore) +
  geom_point(aes(y=m_rsi,x=zscore,colour=Description),
             position = position_jitter(width = .25), size = 4, shape=20) + 
   scale_color_manual(values = c("#125016", "#1b7821", "#24a02c",
                                "#125016", "#000000", "#e7b416",
                                "#CC3232","#a32828","#7a1e1e")) +
  theme_bw() + labs(title="Male",y="Modified Reactive Strength Index")

ggplotly(p6)

p7 <- Df_female_mrsi %>% ggplot(y=m_rsi,x=zscore) +
  geom_point(aes(y=m_rsi,x=zscore,colour=Description),
             position = position_jitter(width = .25), size = 4, shape=20) + 
   scale_color_manual(values = c("#125016", "#1b7821", "#24a02c",
                                "#125016", "#000000", "#e7b416",
                                "#CC3232","#a32828","#7a1e1e")) +
  theme_bw() + labs(title="Female",y="Modified Reactive Strength Index")

ggplotly(p7)

ggarrange(p6,p7) + theme_prism()

ggsave("normative_mrsi_bysex.png")

df_jump_reduced_numerical <- Df[, 10:83]
# Convert all columns to numeric
df_jump_reduced_numerical <- data.frame(lapply(df_jump_reduced_numerical, function(x) as.numeric(as.character(x))))

df_jump_reduced_numerical <- scale(df_jump_reduced_numerical)
corr_matrix <- cor(df_jump_reduced_numerical, use="complete.obs")
ggcorrplot(corr_matrix)

data.pca <- princomp(corr_matrix)
summary(data.pca)

## Importance of components:
##                           Comp.1     Comp.2    Comp.3     Comp.4     Comp.5
## Standard deviation     3.0680060 1.09561279 0.8344215 0.55458296 0.48550661
## Proportion of Variance 0.7698595 0.09817779 0.0569469 0.02515545 0.01927922
## Cumulative Proportion  0.7698595 0.86803730 0.9249842 0.95013965 0.96941886
##                            Comp.6      Comp.7      Comp.8      Comp.9
## Standard deviation     0.35901072 0.291021368 0.253517771 0.168837494
## Proportion of Variance 0.01054178 0.006927058 0.005256733 0.002331508
## Cumulative Proportion  0.97996064 0.986887700 0.992144433 0.994475941
##                            Comp.10     Comp.11      Comp.12     Comp.13
## Standard deviation     0.153957622 0.128868238 0.1061731509 0.086011428
## Proportion of Variance 0.001938659 0.001358285 0.0009219948 0.000605078
## Cumulative Proportion  0.996414600 0.997772885 0.9986948795 0.999299958
##                             Comp.14      Comp.15      Comp.16      Comp.17
## Standard deviation     0.0621456732 0.0502719912 3.049184e-02 2.496177e-02
## Proportion of Variance 0.0003158791 0.0002067051 7.604423e-05 5.096239e-05
## Cumulative Proportion  0.9996158366 0.9998225418 9.998986e-01 9.999495e-01
##                             Comp.18      Comp.19      Comp.20      Comp.21
## Standard deviation     0.0189660860 9.801895e-03 8.839526e-03 4.885771e-03
## Proportion of Variance 0.0000294208 7.858130e-06 6.390826e-06 1.952384e-06
## Cumulative Proportion  0.9999789692 9.999868e-01 9.999932e-01 9.999952e-01
##                             Comp.22      Comp.23      Comp.24      Comp.25
## Standard deviation     4.421451e-03 2.898816e-03 2.737587e-03 2.393806e-03
## Proportion of Variance 1.598927e-06 6.872905e-07 6.129638e-07 4.686805e-07
## Cumulative Proportion  9.999968e-01 9.999975e-01 9.999981e-01 9.999985e-01
##                             Comp.26      Comp.27      Comp.28      Comp.29
## Standard deviation     2.232892e-03 1.649348e-03 1.631149e-03 1.282867e-03
## Proportion of Variance 4.077882e-07 2.224967e-07 2.176137e-07 1.346053e-07
## Cumulative Proportion  9.999989e-01 9.999992e-01 9.999994e-01 9.999995e-01
##                             Comp.30      Comp.31      Comp.32      Comp.33
## Standard deviation     1.233158e-03 1.113051e-03 9.437413e-04 8.585304e-04
## Proportion of Variance 1.243759e-07 1.013280e-07 7.284588e-08 6.028516e-08
## Cumulative Proportion  9.999996e-01 9.999997e-01 9.999998e-01 9.999999e-01
##                             Comp.34      Comp.35      Comp.36      Comp.37
## Standard deviation     7.880441e-04 5.375779e-04 4.296498e-04 3.714300e-04
## Proportion of Variance 5.079256e-08 2.363643e-08 1.509831e-08 1.128374e-08
## Cumulative Proportion  9.999999e-01 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.38      Comp.39      Comp.40      Comp.41
## Standard deviation     2.849653e-04 2.643565e-04 2.096511e-04 1.390099e-04
## Proportion of Variance 6.641759e-09 5.715828e-09 3.594954e-09 1.580487e-09
## Cumulative Proportion  1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.42      Comp.43      Comp.44      Comp.45
## Standard deviation     1.001257e-04 7.351008e-05 5.622760e-05 4.479226e-05
## Proportion of Variance 8.199548e-10 4.419701e-10 2.585819e-10 1.640987e-10
## Cumulative Proportion  1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.46      Comp.47      Comp.48      Comp.49
## Standard deviation     3.544259e-05 2.504252e-05 2.199320e-05 1.882599e-05
## Proportion of Variance 1.027425e-10 5.129265e-11 3.956179e-11 2.898776e-11
## Cumulative Proportion  1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.50      Comp.51      Comp.52      Comp.53
## Standard deviation     1.327168e-05 1.113816e-05 6.596218e-06 6.296993e-06
## Proportion of Variance 1.440625e-11 1.014672e-11 3.558681e-12 3.243138e-12
## Cumulative Proportion  1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.54      Comp.55      Comp.56      Comp.57
## Standard deviation     4.662557e-06 2.913844e-06 2.110931e-06 1.300447e-06
## Proportion of Variance 1.778064e-12 6.944350e-13 3.644579e-13 1.383198e-13
## Cumulative Proportion  1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.58      Comp.59      Comp.60      Comp.61
## Standard deviation     8.896039e-07 7.563176e-07 1.074345e-07 2.348365e-08
## Proportion of Variance 6.472803e-14 4.678509e-14 9.440325e-16 4.510559e-17
## Cumulative Proportion  1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.62      Comp.63      Comp.64      Comp.65
## Standard deviation     1.920058e-08 1.447768e-08 4.014930e-09 3.870902e-09
## Proportion of Variance 3.015281e-17 1.714341e-17 1.318424e-18 1.225528e-18
## Cumulative Proportion  1.000000e+00 1.000000e+00 1.000000e+00 1.000000e+00
##                             Comp.66 Comp.67 Comp.68 Comp.69 Comp.70 Comp.71
## Standard deviation     1.383843e-09       0       0       0       0       0
## Proportion of Variance 1.566292e-19       0       0       0       0       0
## Cumulative Proportion  1.000000e+00       1       1       1       1       1
##                        Comp.72 Comp.73 Comp.74
## Standard deviation           0       0       0
## Proportion of Variance       0       0       0
## Cumulative Proportion        1       1       1

# Get the loadings for the first four principal components
data.pca$loadings[, 3:3]

##                        jump_momentum                countermovement_depth 
##                          0.024332631                         -0.244210680 
##                          braking_rfd                            stiffness 
##                          0.024303734                          0.147711764 
##            force_at_min_displacement   relative_force_at_min_displacement 
##                          0.035015446                          0.031117410 
##                    avg_braking_force           avg_relative_braking_force 
##                         -0.028617452                         -0.062505970 
##                   peak_braking_force          peak_relative_braking_force 
##                          0.035113643                          0.031022853 
##                 avg_propulsive_force        avg_relative_propulsive_force 
##                         -0.066381216                         -0.189357447 
##                peak_propulsive_force       peak_relative_propulsive_force 
##                          0.058370298                          0.063028323 
##                    unweighting_phase            unweighting_phase_percent 
##                          0.217355884                          0.043293166 
##                        braking_phase                braking_phase_percent 
##                          0.206551556                          0.011951446 
##                     propulsive_phase             propulsive_phase_percent 
##                          0.280484119                         -0.082190032 
##                      time_to_takeoff                  braking_net_impulse 
##                          0.279054308                         -0.016667488 
##               propulsive_net_impulse                     positive_impulse 
##                          0.024459341                          0.152586091 
##                 positive_net_impulse                        impulse_ratio 
##                          0.010573646                          0.120819294 
##                 avg_braking_velocity                peak_braking_velocity 
##                          0.003316029                          0.039767062 
##              avg_propulsive_velocity                        peak_velocity 
##                          0.009376006                          0.020263261 
##                    avg_braking_power           avg_relative_braking_power 
##                          0.008652310                          0.020075888 
##                   peak_braking_power          peak_relative_braking_power 
##                          0.004476509                          0.014517207 
##                 avg_propulsive_power        avg_relative_propulsive_power 
##                         -0.037433011                         -0.080679052 
##       peak_relative_propulsive_power                peak_propulsive_power 
##                         -0.047497075                         -0.014039812 
##               l_r_peak_braking_force     left_force_at_peak_braking_force 
##                         -0.192176755                          0.002601037 
##    right_force_at_peak_braking_force                l_r_avg_braking_force 
##                          0.064599522                         -0.163146051 
##               left_avg_braking_force              right_avg_braking_force 
##                         -0.062591812                          0.005625372 
##            l_r_peak_propulsive_force  left_force_at_peak_propulsive_force 
##                         -0.141831830                          0.032414482 
## right_force_at_peak_propulsive_force             l_r_avg_propulsive_force 
##                          0.081170755                         -0.199686477 
##            left_avg_propulsive_force           right_avg_propulsive_force 
##                         -0.092645405                         -0.039740457 
##                  l_r_avg_braking_rfd                 left_avg_braking_rfd 
##                         -0.132578791                          0.003756894 
##                right_avg_braking_rfd            l_r_braking_impulse_index 
##                          0.042169997                         -0.161786987 
##         l_r_propulsive_impulse_index                time_to_stabilization 
##                         -0.199188176                          0.134159710 
##                    landing_stiffness                   peak_landing_force 
##                          0.121752721                         -0.121920208 
##                    avg_landing_force          relative_peak_landing_force 
##                         -0.010515056                         -0.152743952 
##               l_r_peak_landing_force     left_force_at_peak_landing_force 
##                         -0.080023037                         -0.145892271 
##    right_force_at_peak_landing_force                l_r_avg_landing_force 
##                         -0.059586803                         -0.089471715 
##               left_avg_landing_force              right_avg_landing_force 
##                         -0.061192527                          0.036863275 
##            l_r_landing_impulse_index                                m_rsi 
##                         -0.126381572                         -0.161221275 
##                      braking_impulse             relative_braking_impulse 
##                          0.176740309                          0.223412058 
##         relative_braking_net_impulse                   propulsive_impulse 
##                         -0.040848792                          0.122064479 
##          relative_propulsive_impulse      relative_propulsive_net_impulse 
##                          0.251324940                          0.018289322

fviz_eig(data.pca, addlabels = TRUE)

# Graph of the variables
fviz_pca_var(data.pca, col.var = "black")

fviz_cos2(data.pca, choice = "var", axes = 1:2)

fviz_pca_var(data.pca, col.var = "cos2",
            gradient.cols = c("black", "orange", "green"),
            repel = TRUE)

# Assume that 'data.pca' is your PCA result
# Set the threshold
threshold = 0.2

# Get the loadings for the first four principal components
loadings <- data.pca$loadings[, 1:3]

# Find the names of variables where the absolute value of the loading is greater than the threshold in any of the first four components
important_variables <- which(apply(loadings, 1, function(x) any(abs(x) > threshold)))

# If your data.pca was generated by prcomp or princomp, use the following code to extract important variables
# important_variables <- rownames(loadings)[apply(loadings, 1, function(x) any(abs(x) > threshold))]

# Assume that your dataframe is called df
# Create a new dataframe with only these variables
new_df <- Df[, important_variables]

# Print out the new dataframe
print(new_df)

## # A tibble: 115 × 8
##    region relative_force_at_min_…¹ avg_relative_braking…² peak_relative_brakin…³
##    <chr>                     <dbl>                  <dbl>                  <dbl>
##  1 PRK                        248.                   186.                   249.
##  2 SBH                        265.                   178.                   265.
##  3 JHR                        253.                   174.                   253.
##  4 PLS                        220.                   186.                   220.
##  5 MLK                        294.                   192.                   294.
##  6 MLK                        463.                   248.                   464.
##  7 SWK                        240.                   181.                   240.
##  8 SGR                        265.                   186.                   264.
##  9 PHG                        266.                   170.                   265.
## 10 PLS                        254.                   157.                   254.
## # ℹ 105 more rows
## # ℹ abbreviated names: ¹​relative_force_at_min_displacement,
## #   ²​avg_relative_braking_force, ³​peak_relative_braking_force
## # ℹ 4 more variables: avg_relative_propulsive_force <dbl>,
## #   avg_propulsive_velocity <dbl>, l_r_peak_landing_force <dbl>,
## #   l_r_avg_landing_force <dbl>

library(table1)
table1(~ age + jump_height + m_rsi  | sex, data= Df, overall=FALSE)

	Female (N=48)	Male (N=67)
age
Mean (SD)	17.4 (2.76)	18.1 (2.36)
Median [Min, Max]	17.0 [13.0, 23.0]	18.0 [13.0, 23.0]
jump_height
Mean (SD)	0.371 (0.0480)	0.494 (0.0582)
Median [Min, Max]	0.376 [0.268, 0.489]	0.492 [0.358, 0.679]
m_rsi
Mean (SD)	0.437 (0.0862)	0.557 (0.102)
Median [Min, Max]	0.437 [0.261, 0.621]	0.564 [0.245, 0.816]

# Fit the linear mixed model
aov.model <- aov(jump_height ~ sex + age + events, data = Df)

# Print the model summary
summary(aov.model)

##              Df Sum Sq Mean Sq F value Pr(>F)    
## sex           1 0.4201  0.4201 146.864 <2e-16 ***
## age           1 0.0137  0.0137   4.789 0.0308 *  
## events        2 0.0035  0.0018   0.618 0.5411    
## Residuals   110 0.3147  0.0029                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

```